Introduction

Gun violence is a significant issue in the city of Philadelphia. There were 771 shooting incidents in the city in 2024 per the Philadelphia Police Department. The goal of this project is to utilize the data provided by the Philadelphia Police Department via OpenDataPhilly on shooting incidents to build a picture of the fatal shooting incidents occurring in the city. Who is most likely to be the victim of a fatal incident (age, race, sex)? What conditions of an incident make it more likely to be fatal (wound type)? When and where are the most fatal incidents occurring (time of day/year, location)?

Understanding the trends of fatal gun violence incidents can assist in painting a clearer picture of the problem and point towards potential solutions. Understanding the time and location of the most dangerous shootings provides guidance for police on where to focus their efforts. Understanding who is most at risk of a fatal shooting guides social and health services towards populations that may benefit from interventions to assist shooting victims.

Data Overview

The data is collected by the Philadelphia Police Department and provided via OpenDataPhilly. Below, we will load the data and gather the available variables for answering our questions.

Looking at the available fields we have a rich potential set of data to answer our questions. The age, race, sex and latino fields can provide a profile of the most frequent shooting victims. Wound can inform which types of injuries are most often fatal. Location, time and date can be used to determine when and where fatal incidents occur most often. These fields will be compared to the fatal column in order to determine fatality incidence.

# Install necessary packages (run this only once)
# install.packages(c("tidyverse", "sf", "leaflet", "ggplot2", "viridis", "readr", "sp", "leaflet.extras", "readxl", "htmltools", "htmlwidgets", "shiny", "dplyr"))

#install.packages("RColorBrewer")

# Load the packages
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.4     ✔ tidyr     1.3.1
## ✔ purrr     1.0.4     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(sf)
## Linking to GEOS 3.13.0, GDAL 3.10.1, PROJ 9.5.1; sf_use_s2() is TRUE
library(leaflet)
library(ggplot2)
library(viridis)
## Loading required package: viridisLite
library(readr)
library(sp)
library(leaflet.extras)
library(htmltools)
library(htmlwidgets)
library(shiny)
library(dplyr)
library(RColorBrewer)


library(tidyverse)
library(readxl)
shootings.raw = read_delim(
  "./philly_shootings.csv",
  delim=",",
  show_col_types = FALSE,
) %>% glimpse()
## Rows: 16,732
## Columns: 25
## $ the_geom             <chr> "0101000020E61000001049CE6BA7C652C0F0FF5490EA0344…
## $ the_geom_webmercator <chr> "0101000020110F0000790F72E295E45FC1A02823583D9452…
## $ objectid             <dbl> 634246, 634247, 634248, 634249, 634250, 634251, 6…
## $ year                 <dbl> 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2…
## $ dc_key               <dbl> 201502062409, 201503069342, 201512041905, 2015120…
## $ code                 <dbl> 3400, 400, 400, 400, 400, 400, 400, 400, 3400, 40…
## $ date_                <dttm> 2015-09-26, 2015-11-08, 2015-05-30, 2015-12-11, …
## $ time                 <time> 05:07:00, 12:50:00, 14:11:00, 18:00:00, 19:40:00…
## $ race                 <chr> "B", "A", "B", "B", "B", "B", "B", "B", "B", "B",…
## $ sex                  <chr> "M", "M", "M", "M", "F", "F", "F", "M", "F", "M",…
## $ age                  <dbl> 43, 36, 27, 23, 15, 21, 26, 33, 26, 43, 21, 35, 2…
## $ wound                <chr> "Chest", "Arm", "Stomach", "Leg", "Leg", "Leg", "…
## $ officer_involved     <chr> "N", "N", "N", "N", "N", "N", "N", "N", "N", "N",…
## $ offender_injured     <chr> "N", "N", "N", "N", "N", "N", "N", "N", "N", "N",…
## $ offender_deceased    <chr> "N", "N", "N", "N", "N", "N", "N", "N", "N", "N",…
## $ location             <chr> "700 BLOCK Adams Ave", "1700 BLOCK S 5th St", "69…
## $ latino               <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ point_x              <dbl> -75.10397, -75.15405, -75.23670, -75.22840, -75.1…
## $ point_y              <dbl> 40.03060, 39.92726, 39.91765, 39.91845, 40.04942,…
## $ dist                 <dbl> 2, 3, 12, 12, 14, 14, 14, 14, 14, 14, 14, 14, 14,…
## $ inside               <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ outside              <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ fatal                <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ lat                  <dbl> 40.03060, 39.92726, 39.91765, 39.91845, 40.04942,…
## $ lng                  <dbl> -75.10397, -75.15405, -75.23670, -75.22840, -75.1…

Exploratory Data Analysis

Handling NAs

First, we look at NAs to see where we may have missing data that needs to be cleaned up. Looking at the summary of NAs a group of 153 NAs shows up across multiple columns including fatal. We will investigate those to clean them up since they will prevent us from accurately summing the fatalities.

shootings.raw %>% summarise_all(~ sum(is.na(.)))

Based on this, it looks like certain columns, including fatal, are left blank in any incident that involves an officer. In these cases the only fatality data provided is offender deceased. However, this would leave out any other fatalities that occurred during the incident since not everyone injured in a shooting incident is necessarily an offender.

Critical data for our analysis is missing from all of these officer involved shootings - time, race, age, wound, latino and, most importantly, fatal. There are 40 fatal incidents with an officer involved for which there is incomplete data.

shootings.raw %>% filter(if_any(fatal, is.na))
shootings.raw %>% filter(officer_involved == 'Y')
shootings.officer = shootings.raw %>% filter(officer_involved == 'Y')

shootings.officer %>% filter(offender_deceased == 'Y')

We can make an adjusted table, adding a total_fatality column that includes offender_deceased. It will not be usable for the questions of the profile of the most frequent shooting victims and their wound types, or the time of day of incidents, since that data is missing, but it can be used for the geolocation and date questions.

This introduces some potential bias into our data because it is possible that we are excluding non-offender fatalities. It is unknown to us if those are excluded or just did not occur. However, it is more complete than not including these 40 fatalities at all.

This gives us a column, total_fatality, with no NAs that can be summed to a total of 3494 fatalities between 2015 and 2025.

shootings.adjusted = shootings.raw

shootings.adjusted$offender_fatal <- ifelse(shootings.raw$offender_deceased == "Y", 1, 0)

shootings.adjusted <- shootings.adjusted %>% 
                      mutate(fatal = replace_na(fatal, 0))

shootings.adjusted$total_fatality <- shootings.adjusted$offender_fatal + shootings.adjusted$fatal

shootings.adjusted
sum(shootings.adjusted$total_fatality)
## [1] 3496

Examining Data by Year

Looking at fatality totals by year we can see that partial data is included for 2025. This will obviously be skewed on an annual basis since 2025 is not complete yet. We can see here that fatal shooting incidents rose from 2015 to 2021, peaked in 2021, then fell from 2021-2024.

fatal.byyear = aggregate(shootings.adjusted['total_fatality'], by=shootings.adjusted['year'], sum)

fatal.15to24 <- fatal.byyear %>%
  filter(year != 2025)

ggplot(data = fatal.15to24, aes(x = year, y = total_fatality)) +
  ggtitle("Total Fatal Shooting Victims by Year 2015-2024") +
  geom_col() +
  scale_x_continuous(breaks = 2015:2024)

AK - Geospatial Analysis

Clean and prepare the data

# View structure
glimpse(shootings.adjusted)
## Rows: 16,732
## Columns: 27
## $ the_geom             <chr> "0101000020E61000001049CE6BA7C652C0F0FF5490EA0344…
## $ the_geom_webmercator <chr> "0101000020110F0000790F72E295E45FC1A02823583D9452…
## $ objectid             <dbl> 634246, 634247, 634248, 634249, 634250, 634251, 6…
## $ year                 <dbl> 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2…
## $ dc_key               <dbl> 201502062409, 201503069342, 201512041905, 2015120…
## $ code                 <dbl> 3400, 400, 400, 400, 400, 400, 400, 400, 3400, 40…
## $ date_                <dttm> 2015-09-26, 2015-11-08, 2015-05-30, 2015-12-11, …
## $ time                 <time> 05:07:00, 12:50:00, 14:11:00, 18:00:00, 19:40:00…
## $ race                 <chr> "B", "A", "B", "B", "B", "B", "B", "B", "B", "B",…
## $ sex                  <chr> "M", "M", "M", "M", "F", "F", "F", "M", "F", "M",…
## $ age                  <dbl> 43, 36, 27, 23, 15, 21, 26, 33, 26, 43, 21, 35, 2…
## $ wound                <chr> "Chest", "Arm", "Stomach", "Leg", "Leg", "Leg", "…
## $ officer_involved     <chr> "N", "N", "N", "N", "N", "N", "N", "N", "N", "N",…
## $ offender_injured     <chr> "N", "N", "N", "N", "N", "N", "N", "N", "N", "N",…
## $ offender_deceased    <chr> "N", "N", "N", "N", "N", "N", "N", "N", "N", "N",…
## $ location             <chr> "700 BLOCK Adams Ave", "1700 BLOCK S 5th St", "69…
## $ latino               <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ point_x              <dbl> -75.10397, -75.15405, -75.23670, -75.22840, -75.1…
## $ point_y              <dbl> 40.03060, 39.92726, 39.91765, 39.91845, 40.04942,…
## $ dist                 <dbl> 2, 3, 12, 12, 14, 14, 14, 14, 14, 14, 14, 14, 14,…
## $ inside               <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ outside              <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ fatal                <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ lat                  <dbl> 40.03060, 39.92726, 39.91765, 39.91845, 40.04942,…
## $ lng                  <dbl> -75.10397, -75.15405, -75.23670, -75.22840, -75.1…
## $ offender_fatal       <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
## $ total_fatality       <dbl> 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
summary(shootings.adjusted)
##    the_geom         the_geom_webmercator    objectid           year     
##  Length:16732       Length:16732         Min.   :634246   Min.   :2015  
##  Class :character   Class :character     1st Qu.:638429   1st Qu.:2018  
##  Mode  :character   Mode  :character     Median :642612   Median :2020  
##                                          Mean   :642612   Mean   :2020  
##                                          3rd Qu.:646794   3rd Qu.:2022  
##                                          Max.   :650977   Max.   :2025  
##                                                                         
##      dc_key               code            date_                       
##  Min.   :1.502e+03   Min.   : 100.0   Min.   :2015-01-01 00:00:00.00  
##  1st Qu.:2.018e+11   1st Qu.: 300.0   1st Qu.:2018-04-02 00:00:00.00  
##  Median :2.020e+11   Median : 400.0   Median :2020-09-25 00:00:00.00  
##  Mean   :2.001e+11   Mean   : 428.5   Mean   :2020-05-14 23:53:22.39  
##  3rd Qu.:2.022e+11   3rd Qu.: 400.0   3rd Qu.:2022-06-27 06:00:00.00  
##  Max.   :2.025e+11   Max.   :3700.0   Max.   :2025-04-28 00:00:00.00  
##                      NA's   :167                                      
##      time              race               sex                 age        
##  Length:16732      Length:16732       Length:16732       Min.   :  0.00  
##  Class1:hms        Class :character   Class :character   1st Qu.: 21.00  
##  Class2:difftime   Mode  :character   Mode  :character   Median : 27.00  
##  Mode  :numeric                                          Mean   : 29.34  
##                                                          3rd Qu.: 35.00  
##                                                          Max.   :117.00  
##                                                          NA's   :274     
##     wound           officer_involved   offender_injured   offender_deceased 
##  Length:16732       Length:16732       Length:16732       Length:16732      
##  Class :character   Class :character   Class :character   Class :character  
##  Mode  :character   Mode  :character   Mode  :character   Mode  :character  
##                                                                             
##                                                                             
##                                                                             
##                                                                             
##    location             latino         point_x          point_y     
##  Length:16732       Min.   :0.000   Min.   :-75.27   Min.   :39.88  
##  Class :character   1st Qu.:0.000   1st Qu.:-75.19   1st Qu.:39.97  
##  Mode  :character   Median :0.000   Median :-75.16   Median :39.99  
##                     Mean   :0.123   Mean   :-75.16   Mean   :39.99  
##                     3rd Qu.:0.000   3rd Qu.:-75.13   3rd Qu.:40.02  
##                     Max.   :1.000   Max.   :-74.96   Max.   :40.13  
##                     NA's   :153     NA's   :33       NA's   :33     
##       dist           inside           outside           fatal       
##  Min.   : 1.00   Min.   :0.00000   Min.   :0.0000   Min.   :0.0000  
##  1st Qu.:15.00   1st Qu.:0.00000   1st Qu.:1.0000   1st Qu.:0.0000  
##  Median :22.00   Median :0.00000   Median :1.0000   Median :0.0000  
##  Mean   :20.91   Mean   :0.06333   Mean   :0.9367   Mean   :0.2066  
##  3rd Qu.:25.00   3rd Qu.:0.00000   3rd Qu.:1.0000   3rd Qu.:0.0000  
##  Max.   :77.00   Max.   :1.00000   Max.   :1.0000   Max.   :1.0000  
##  NA's   :2       NA's   :153       NA's   :153                      
##       lat             lng         offender_fatal     total_fatality  
##  Min.   :39.88   Min.   :-75.27   Min.   :0.000000   Min.   :0.0000  
##  1st Qu.:39.97   1st Qu.:-75.19   1st Qu.:0.000000   1st Qu.:0.0000  
##  Median :39.99   Median :-75.16   Median :0.000000   Median :0.0000  
##  Mean   :39.99   Mean   :-75.16   Mean   :0.002391   Mean   :0.2089  
##  3rd Qu.:40.02   3rd Qu.:-75.13   3rd Qu.:0.000000   3rd Qu.:0.0000  
##  Max.   :40.13   Max.   :-74.96   Max.   :1.000000   Max.   :1.0000  
##  NA's   :28      NA's   :28
shootings.spatial  <- shootings.adjusted %>%
  filter(!is.na(lat) & !is.na(lng)) %>%
  rename(latitude = lat, longitude = lng)

Convert to spatial data

shootings_sf <- st_as_sf(shootings.spatial , coords = c("longitude", "latitude"), crs = 4326)

Plotting Static Map

ggplot(data = shootings_sf) +
  geom_sf(alpha = 0.3, color = "red") +
  theme_minimal() +
  labs(title = "Philadelphia Shooting Incidents",
       subtitle = "Based on geospatial coordinates",
       caption = "Source: philly_shootings.csv")

Interactive Map Using Leaflet

# Load necessary libraries
# Ensure you have these: dplyr, sf, leaflet, RColorBrewer, grDevices, htmltools, htmlwidgets, jsonlite
library(dplyr)
library(sf)
library(leaflet)
library(RColorBrewer)
library(grDevices)
library(htmltools)
library(htmlwidgets)
library(jsonlite)

# Ensure shootings_sf is an sf object and has relevant columns
# This assumes 'shootings_sf' is already loaded and processed.
# Replace the placeholder below with your actual data loading and preparation.
if (!exists("shootings_sf") || !"sf" %in% class(shootings_sf)) {
 message("Warning: 'shootings_sf' not found or not an sf object. Using dummy data for example structure.")
 # Dummy data for example structure:
 shootings_sf <- data.frame(
    year = rep(2020:2021, each = 2),
    date_ = as.Date(c("2020-01-01", "2020-02-01", "2021-01-15", "2021-03-10")),
    location = paste("Location", 1:4),
    total_fatality = c(0, 1, 0, 1),
    lat = runif(4, 39.9, 40.1),
    lng = runif(4, -75.2, -75.1)
 ) %>%
 st_as_sf(coords = c("lng", "lat"), crs = 4326)
}

if (!"total_fatality" %in% names(shootings_sf) || !"year" %in% names(shootings_sf) || !"date_" %in% names(shootings_sf) || !"location" %in% names(shootings_sf) ) {
  stop("shootings_sf is missing one or more required columns: year, date_, location, total_fatality.")
}

if (!"unique_id" %in% names(shootings_sf)) {
  shootings_sf <- shootings_sf %>%
    dplyr::mutate(unique_id = paste0("marker-", dplyr::row_number(), "-", year))
}

# Create unique list of years
years <- if ("year" %in% names(shootings_sf)) sort(unique(shootings_sf$year)) else character(0)

# Create color palette
if (length(years) > 0) {
  if (length(years) <= 8) {
    year_palette_colors <- RColorBrewer::brewer.pal(max(3, length(years)), "Dark2")
    if (length(years) < 3) {
      year_palette_colors <- year_palette_colors[1:length(years)]
    }
  } else {
    year_palette_colors <- grDevices::colorRampPalette(RColorBrewer::brewer.pal(min(length(years), 9), "Dark2"))(length(years))
  }
  year_colors <- leaflet::colorFactor(palette = year_palette_colors, domain = years)
} else {
  year_colors <- NULL
  warning("Year column not found or no unique years in shootings_sf. Map may not display year layers correctly.")
}

# Create base map centered on Philadelphia
center_lon <- -75.1652
center_lat <- 39.9526

m <- leaflet::leaflet(shootings_sf) %>%
  leaflet::setView(lng = center_lon, lat = center_lat, zoom = 11) %>%
  leaflet::addProviderTiles(leaflet::providers$CartoDB.Positron) %>%
  htmlwidgets::onRender("function(el, x) { el.classList.add('safegraph-map'); }")

# Add circle markers by year
if (length(years) > 0 && !is.null(year_colors)) {
  for (yr in years) {
    year_data <- shootings_sf %>% dplyr::filter(year == yr)
    if (nrow(year_data) > 0) {
      m <- m %>% leaflet::addCircleMarkers(
        data = year_data,
        radius = 5,
        color = ~year_colors(year),
        fillOpacity = 0.7,
        stroke = FALSE,
        group = as.character(yr),
        popup = ~paste(
          "<b>Year:</b>", htmltools::htmlEscape(year), "<br>",
          "<b>Date:</b>", htmltools::htmlEscape(date_), "<br>",
          "<b>Location:</b>", htmltools::htmlEscape(location), "<br>",
          "<b>Fatality:</b>", ifelse(total_fatality > 0, "Yes", "No"),
          "<b>Wound:</b>", htmltools::htmlEscape(wound)
        ),
        layerId = ~unique_id
      )
    }
  }
}

# Add layer control
if (length(years) > 0) {
  m <- m %>% leaflet::addLayersControl(
    overlayGroups = as.character(years),
    options = leaflet::layersControlOptions(collapsed = FALSE)
  )
}

# JavaScript for a single "Select/Deselect All" toggle button using paste0
year_groups_json_str <- jsonlite::toJSON(as.character(years)) # Prepare JSON string once

select_toggle_js <- paste0(
  "function(el, x, data) {\n",
  "  setTimeout(function() {\n",
  "    var yearGroups = ", year_groups_json_str, ";\n",
  "    if (!yearGroups || yearGroups.length === 0) { console.warn('No year groups for toggle button.'); return; }\n",
  "    var layersControlContainer = el.querySelector('.leaflet-control-layers-overlays');\n",
  "    if (!layersControlContainer) { console.warn('Layers control container not found for toggle button.'); return; }\n",
  "\n",
  "    var toggleButton = L.DomUtil.create('button', 'leaflet-control-layers-selector safegraph-button');\n",
  "    toggleButton.style.marginBottom = '10px';\n",
  "\n",
  "    function getCurrentAllCheckedState() {\n",
  "      if (yearGroups.length === 0) return false;\n",
  "      for (var i = 0; i < yearGroups.length; i++) {\n",
  "        var groupStr = String(yearGroups[i]);\n",
  "        var checkbox = layersControlContainer.querySelector('input[type=\"checkbox\"][value=\"' + groupStr + '\"]');\n",
  "        if (!checkbox || !checkbox.checked) return false;\n",
  "      }\n",
  "      return true;\n",
  "    }\n",
  "\n",
  "    function updateButtonText() {\n",
  "      if (getCurrentAllCheckedState()) {\n",
  "        toggleButton.innerHTML = 'Deselect All Years';\n",
  "      } else {\n",
  "        toggleButton.innerHTML = 'Select All Years';\n",
  "      }\n",
  "    }\n",
  "\n",
  "    toggleButton.onclick = function() {\n",
  "      var currentlyAllChecked = getCurrentAllCheckedState();\n",
  "      var targetState = !currentlyAllChecked;\n",
  "      yearGroups.forEach(function(group) {\n",
  "        var groupStr = String(group);\n",
  "        var checkbox = layersControlContainer.querySelector('input[type=\"checkbox\"][value=\"' + groupStr + '\"]');\n",
  "        if (checkbox) {\n",
  "          if (checkbox.checked !== targetState) {\n",
  "            checkbox.click();\n",
  "          }\n",
  "        }\n",
  "      });\n",
  "      updateButtonText();\n",
  "    };\n",
  "\n",
  "    layersControlContainer.parentNode.insertBefore(toggleButton, layersControlContainer);\n",
  "    updateButtonText(); /* Set initial button text */\n",
  "  }, 250);\n",
  "}\n"
)

m <- m %>% htmlwidgets::onRender(select_toggle_js)

# Add legend
if (length(years) > 0 && !is.null(year_colors)) {
  m <- m %>% leaflet::addLegend(
    position = "bottomleft",
    pal = year_colors,
    values = shootings_sf$year,
    title = "Shooting Year",
    opacity = 1
  )
}

# CSS Styling using paste0
css_styling_js <- paste0(
  "function(el, x) {\n",
  "  var css = `\n", # Start of JS template literal
  "    .safegraph-map { border-radius: 8px; box-shadow: 0 2px 8px rgba(0,0,0,0.15); overflow: hidden; }\n",
  "    .safegraph-button { background-color: #f0f0f0; border: 1px solid #ccc; border-radius: 4px; padding: 4px 8px; font-size: 12px; cursor: pointer; transition: background-color 0.3s ease; font-family: 'Arial', sans-serif; color: #333; display: block; width: 100%; text-align: center; white-space: nowrap; box-sizing: border-box; }\n",
  "    .safegraph-button:hover { background-color: #e0e0e0; }\n",
  "    .leaflet-control-layers-list input[type='checkbox'] { margin-right: 4px; }\n",
  "    .leaflet-popup-content { font-size: 12px; line-height: 1.4; }\n",
  "    .leaflet-popup-content b { font-weight: 600; color: #2c3e50; }\n",
  "    .leaflet-control-layers { background-color: rgba(255,255,255,0.9); border-radius: 8px; box-shadow: 0 1px 4px rgba(0,0,0,0.1); padding: 6px 10px; font-size: 12px; color: #34495e; }\n",
  "    .leaflet-control-layers-title { font-weight: bold; margin-bottom: 6px; color: #2c3e50; }\n",
  "  `;\n", # End of JS template literal
  "  var style = document.createElement('style');\n",
  "  style.type = 'text/css';\n",
  "  style.appendChild(document.createTextNode(css));\n",
  "  document.head.appendChild(style);\n",
  "  if (!el.classList.contains('safegraph-map')) { el.classList.add('safegraph-map'); }\n",
  "}\n"
)

m <- m %>% htmlwidgets::onRender(css_styling_js)

# Display map
m
# R chunk: Defines map_enhanced_heatmaps (Corrected)

# Define map center
center_lon <- -75.1652
center_lat <- 39.9526

# Prepare base data for the map
# This assumes 'shootings_sf' contains 'officer_involved' and 'offender_deceased' from your raw data processing
shootings_sf_for_map <- shootings_sf %>%
  dplyr::mutate(
    incident_color = ifelse(total_fatality > 0, "red", "orange"),
    fatality_label = ifelse(total_fatality > 0, "Fatal", "Non-Fatal")
  )

# Check for necessary columns for filtering new heatmaps
if (!all(c("officer_involved", "offender_deceased", "total_fatality") %in% names(shootings_sf_for_map))) {
  stop("One or more required columns ('officer_involved', 'offender_deceased', 'total_fatality') are missing from shootings_sf_for_map.")
}
if (!"wound" %in% names(shootings_sf_for_map)) {
  warning("The column 'wound' was not found in the data. Wound information will be missing from popups.")
}


# Create filtered datasets for new heatmap layers
officer_involved_fatal_sf <- shootings_sf_for_map %>%
  dplyr::filter(officer_involved == 'Y' & total_fatality > 0)

officer_involved_non_fatal_sf <- shootings_sf_for_map %>%
  dplyr::filter(officer_involved == 'Y' & total_fatality == 0)

offender_deceased_sf <- shootings_sf_for_map %>%
  dplyr::filter(offender_deceased == 'Y') 

all_officer_involved_sf <- shootings_sf_for_map %>%
  dplyr::filter(officer_involved == 'Y')
  
# Create the Leaflet map
map_enhanced_heatmaps <- leaflet::leaflet(data = shootings_sf_for_map) %>% # Base data
  leaflet::setView(lng = center_lon, lat = center_lat, zoom = 11) %>%
  leaflet::addProviderTiles(leaflet::providers$CartoDB.Positron, group = "Base Map") %>%
  
  # Clustered Markers (as before)
  leaflet::addAwesomeMarkers(
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2], # <--- THIS LINE WAS CORRECTED
    icon = ~leaflet::awesomeIcons(
      icon = 'info-circle', 
      library = 'glyphicon', 
      markerColor = incident_color, 
      iconColor = 'white'
    ),
    popup = ~paste(
      "<b>Year:</b>", htmltools::htmlEscape(year), "<br>",
      "<b>Date:</b>", htmltools::htmlEscape(date_), "<br>",
      "<b>Location:</b>", htmltools::htmlEscape(location), "<br>",
      "<b>Fatality:</b>", ifelse(total_fatality > 0, "Yes", "No"), "<br>",
      "<b>Wound:</b>", if ("wound" %in% names(.)) htmltools::htmlEscape(wound) else "N/A"
    ),
    clusterOptions = leaflet::markerClusterOptions(),
    group = "Shooting Incidents (Clusters)"
  ) %>%
  
  # Heatmap
  leaflet.extras::addHeatmap(
    # Data from leaflet() call: shootings_sf_for_map
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2],
    blur = 20, max = 0.05, radius = 15,
    group = "Heatmap (All Incidents)"
  ) %>%

  # Heatmap for All Fatal Incidents (existing - uses total_fatality)
  leaflet.extras::addHeatmap(
    data = dplyr::filter(shootings_sf_for_map, total_fatality > 0), 
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2],
    blur = 20, max = 0.05, radius = 15,
    group = "Heatmap (All Fatal Incidents)" 
  ) %>%
  
  # Heatmap for Officer Involved & Fatal Incidents
  leaflet.extras::addHeatmap(
    data = officer_involved_fatal_sf,
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2],
    blur = 20, max = 0.05, radius = 15,
    group = "Heatmap (Officer Involved & Fatal)"
  ) %>%
  
  leaflet.extras::addHeatmap(
    data = officer_involved_non_fatal_sf,
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2],
    blur = 20, max = 0.05, radius = 15, 
    group = "Heatmap (Officer Involved & Non-Fatal)"
  ) %>%
  
  # Heatmap for Offender Deceased Incidents
  leaflet.extras::addHeatmap(
    data = offender_deceased_sf,
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2],
    blur = 20, max = 0.05, radius = 15,
    group = "Heatmap (Offender Deceased)"
  ) %>%

  # Heatmap for All Officer Involved Incidents
  leaflet.extras::addHeatmap(
    data = all_officer_involved_sf,
    lng = ~sf::st_coordinates(geometry)[,1], 
    lat = ~sf::st_coordinates(geometry)[,2],
    blur = 20, max = 0.05, radius = 15,
    group = "Heatmap (All Officer Involved)"
  ) %>%
  
  # Legend for Clustered Markers
  leaflet::addLegend(
    position = "bottomright", 
    colors = c("red", "orange"),
    labels = c("Fatal Incident", "Non-Fatal Incident"),
    title = "Incident Type (Markers)",
    opacity = 1,
    group = "Shooting Incidents (Clusters)" 
  ) %>%
  
  # Layers Control to switch between layers
  leaflet::addLayersControl(
    baseGroups = "Base Map",
    overlayGroups = c(
      "Shooting Incidents (Clusters)", 
      "Heatmap (All Incidents)", 
      "Heatmap (All Fatal Incidents)",
      "Heatmap (All Officer Involved)",
      "Heatmap (Officer Involved & Fatal)",
      "Heatmap (Officer Involved & Non-Fatal)", 
      "Heatmap (Offender Deceased)"
    ),
    options = leaflet::layersControlOptions(collapsed = FALSE) 
  ) %>%
  # Hide all heatmap layers by default so clusters are visible first
  leaflet::hideGroup(c(
      "Heatmap (All Incidents)", 
      "Heatmap (All Fatal Incidents)",
      "Heatmap (All Officer Involved)",
      "Heatmap (Officer Involved & Fatal)",
      "Heatmap (Officer Involved & Non-Fatal)",
      "Heatmap (Offender Deceased)"
  ))

# Display the map
map_enhanced_heatmaps